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Abstract 

We present an asymptotic analysis of the Gunn effect in a drift- 
diffusion model — including electric-field-dependent generation-recombination 
processes — for long samples of strongly compensated p-type Ge at low tem- 
perature and under dc voltage bias. During each Gunn oscillation, there are 
different stages corresponding to the generation, motion and annihilation of 
solitary waves. Each stage may be described by one evolution equation for 
only one degree of freedom (the current density), except for the generation 
of each new wave. The wave generation is a faster process that may be de- 
scribed by solving a semiinfinite canonical problem. As a result of our study 
we have found that (depending on the boundary condition) one or several 
solitary waves may be shed during each period of the oscillation. Examples 
of numerical simulations validating our analysis are included. 
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I. INTRODUCTION 



In recent years a great variety of oscillatory behaviors have been observed in semiconduc- 
tors displaying nonlinear electrical conduction. A particularly interesting system to study 
spatiotemporal phenomena in the laboratory is ultrapure extrinsic cooled bulk p-type Ge 
under voltage bias conditions, . Observed phenomena include time-periodic oscillations 
of the current on a purely resistive external circuit under dc voltage bias due to the periodic 
creation of a solitary wave at the injecting contact, its motion inside the semiconductor and 
its annihilation at the receiving contact 0. Essentially this is the same as the usual Gunn 
effect in n-GaAs [|], except that the solitary wave in p-Ge moves much more slowly than the 
carriers (the case in the usual Gunn effect) due to the generation-recombination dynamics 
of ionized traps which dominate the transport properties. Other experimental observations 
include intermittency near the onset of the oscillatory instability and spatiotemporal 
chaos under dc-fac voltage bias [^. These phenomena have been studied theoretically by 
means of a drift- diffusion model which includes impurity trapping of mobile holes and im- 
pact ionization of neutral acceptors [^|[. In dimensionless form, the model equations are 

8 A r 

-{a-l-A) + P[{a-l)K{E) 



dt p 



- R{E) - [K{E) + R{E)] A] , (1) 

,3-^j.V(E)P^S-, (2) 

1 

E{x,t) dx = (p. (4) 



L Jo 

The unknowns in these equations are the electric field E{x, t), the hole concentration P{x, t), 
the concentration of the ionized acceptors A{x,t) and the total current density J{t). (p is 
the constant dc bias. To solve this system of equations, we need two boundary conditions 
(representing the effect of the contacts at x = and x = L) and two initial conditions. In 
Section |T| we shall prove that we have a well-posed problem with positive solutions for all 
time for the following Dirichlet boundary conditions: 

P(0,t)=Po, P{L,t)=P,. (5) 
(the Pj's are positive numbers). These conditions have been considered before in [|l^. A 



discussion of imperfect boundary conditions for the Gunn effect introducing the concept of 
contact characteristic was first published by H. Kroemer [^. Ohmic boundary conditions in 
which the field at the contact is proportional to the current density may give rise to negative 
current densities and numerical instabilities as time elapses. Eq. (|l|) is a rate equation for 
the ionized acceptors including the processes of photogeneration, impact ionization and 
recombination (coefficients F, K{E) and R{E), respectively). Eq. (Q) is a form of Ampere's 
law saying that the displacement current is equal to the total current density minus the 
drift-diffusion current due to the free holes. (0) is the Poisson equation and Eq. (^ is the 
voltage bias condition. The transport coefficients V{E), K{E) and R{E) depend nonlinearly 
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on the electric field as shown in Fig. |l] §]. The qualitative nature of much of the predicted 
behavior — e.g., the instability of the stationary electric-field profile and the stability of 
propagating high-field domains — depends only on the presence of a negative slope of the 
homogeneous stationary current density j{E) (defined below) over an interval of positive 
fields, P,p!2| and not on the exact form of the underlying coefficients. Other authors have 



stated different forms of the coefficient curves for various reasons. JTSHII Recently, Monte 
Carlo simulation has been used to provide a more rigorous determination of the coefficients 



for a related p-type system, |T6[ however such simulations do not yet exist for the closely 
compensated Ge samples under consideration in this paper. The homogeneous stationary 
current density j{E) is given by the formula 

P, p!7| , p!8| . If the compensation ratio a (the ratio of the acceptor concentration to the donor 
concentration) is only slightly larger than 1 (see [jl9| for the details), j{E) is N-shaped for 
large enough positive fields. Then there is an interval {Em, Em) between the abcissas of the 
maximum [J^Em) = Jm] and the minimum [j(-Em) = jm > 0] oi i{E) for which dj/dE < 
and j{E) > 0. See Fig. |. 

The dimensionless parameters j3, 5, and F//? are very small for the p-Ge and can be set 
to zero in the leading order approximation Then P = J/V{E) from (^, and (|Tp and (^ 
become: 

dE J 

Using Eq. (|^), we may eliminate A in favor of E and then insert the result into (^. We 
obtain the following reduced equation [pO|,p|,pTl : 

d'^E K{E) + R{E) / V'{E) dE 
V{Ey [k{E) + R{E) 'dt 

The reduced model equations (^) and (H) for E{x,t) and J(t) are supplemented with 
the boundary condition at the injecting contact J(t) = PoV[E{0,t)]. According to this 
equation, the electric field at the contact is a nonlinear function of the current. In |^T| it 
was shown that the resulting stable current oscillations are qualitatively like the solutions 
of the same reduced problem with the linear relation: 

E{0,t) = pJ{t), (10) 

(Ohm's law; p > is the dimensionless contact resistivity) and an initial condition E{x,0). 



Although much work has been done on the reduced model problem (^, (^ and (|TOD [see 



T8| and references therein; equivalently we may analyze (^, (|^), @ and (0)], important 
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basic questions concerning the asymptotic description of their stable solutions are still open. 
In this paper we present an asymptotic description of the Gunn effect in p-Ge and compare 
the results with numerical simulations. 

The paper is organized as follows. Section |I| contains a proof of the global existence of 
positive solutions to the full drift-diffusion model with particular boundary conditions and a 
brief review of the known results for the reduced model equations. In Sec. |T|, we introduce 
our asymptotic scaling, derive the outer and inner solutions, and put the pieces together 
thereby describing one period of the solution through the succession of its different stages. 
The delicate problem of shedding new waves through an instability of the boundary layer 
near the injecting contact is discussed in Section Sec. |V| contains our conclusions and 
some open problems. Appendix A contains technical details on the shedding problem. 



II. MODEL SYSTEM AND GLOBAL EXISTENCE AND UNIQUENESS OF 

SOLUTIONS 

We shall next show that, under rather mild assumptions on coefficients and data, the 
system ([|)-(||) has a unique global solution when supplemented with initial conditions 

A{x,0) = Ao{x), P(x,0) = P(x). (11) 

Specifically, we will require functions Aq, P, Pq and Pi (Pq ^-nd Pi may be known functions 
of time, not just constants) to be bounded, whereas functions K{E), R{E) and V{E) will be 
assumed to be C^. Furthermore, we shall impose the following growth condition on K{E): 

\K{E)\ < C{\E\^ + 1) forsome7>0 (12) 

(which is perfectly reasonable in view of the physics behind the model; see [rB,TD|). 

Without loss of generality, we may assume Pq = Pi = in (H), and F = in (|lD. As a 
matter of fact, the term thus dropped in (|l|) is a linear one and therefore harmless from the 
point of view of global existence. For the sake of notational simplicity, we shall also write: 

f{A, E) = {a- l)K{E) - R{E) - [K{E) + R{E)] A. (13) 

Differentiating (^ with respect to the space coordinate x and using (^, we find the equation 
of charge continuity: 

f3Pt = 6P,, - {V{E) P)^ + /3 A, (14) 

which is equivalent to 

[3Pt = 6P,, - {V{E) P)^ + [3 f{A, E) P, (15) 

once (|T]) is substituted into (|np. Let us discuss first local existence of solutions. To this 
end, we observe that on setting ip = — Jq E{C,,t)dC, (the electric potential), equations (^ 
and (|^) read: 

-ip.ee = P - A; ij{0,t) = 0, ij{L,t) = -0L 
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By classical results for Poisson's equation, we then deduce that: 

Lemma 1. If P and A are bounded, then {E^l = lipxxl < C for any t > 0, where C 
depends only on the hounds for P and A. 

On the other hand, integrating in (|^) gives: 

E = r [p(e,t) - A(e,t)] di + Ci, 

JO 



where Ci = Ci{(j), P, A). Hence, equation (15) may be rewritten in the form: 

(3Pt = 6P,, - {V{E) P)^ + /3 f{A, E{P, A)) P (16) 
= 6Pxx + MiP,Px,E,E,,A), 

where the operator Ai is locally Lipschitz in its arguments. To proceed further, we remark 
that an integration in time of (|I|) yields 

A{x, t) = A{x, 0) + /* f[A{x, s), E{P, A)] P{x, s)ds (17) 
Jo 

In view of ([T6| ) and (0), local existence follows now by means of a contractive fixed point 
argument for the following operator: 



T 



/ P \ / p{t)P + \ fl S{t - s)M{P, P., E, E,, A)ds 
^a) U(^,0) + J^fiAix,s),EiP,A)) 



(18) 



where S{t) denotes the heat semigroup in [0, L] with homogeneous Dirichlet conditions. 
The operator T is to be considered as acting on the space of functions (P, A) such that 
P e L°°(0, t*; H^{0, L)) and A G L°°(0, t*- L°°(0, L)), where t* is a sufficiently small positive 
number. On using the standard estimate 

||5(t)P||Hi(0,L) < Ct-^/^\\P\\L2i^,^L) 

for some C > (see for instance I^Sj), we readily derive the contractivity of the operator 
T for < t* ^ 1. This gives at once local existence and uniqueness. We next set out to 
describe our global existence argument. To that end, we shall first prove that the following 
result holds for (|1^): 

Lemma 2. Assume that A{x,t) and P[x) are bounded. Then for any ai > 1 there exists 
C = C{ai) > such that 

^* \P{x,t)\''' dxdt^ 

< C (p|U(t+l) + ||P||oo), (19) 
for any t > for which a solution of our problem exists. 
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We remark in passing that we are denoting by ||v4||oo a bound on the supremum of A. 
To prove Lemma 2, we argue as follows. Set 

u = f P{x, s) ds (20) 

J 

Then u satisfies 

^ du ^ d'^u ^ , , , 



dx 



j\{E{x,s))P{x,s)ds^ . (21) 



If the last terms on the right of ( [TE| ) were absent, standard parabolic theory would yield at 
once 

\utrdxdtj'''' < C (PIU + ||P||oc) (22) 

for any ai > 1, 

where Q = (0, L) x (0, t*) and C = C{L, t*) > 0. That would yield (||) except for the linear 
dependence on t* stated there. To prove Lemma 2 in the general case, we shall assume for 
simplicity that P = 0, and note that a solution of (^T]) with homogeneous initial and side 
conditions may be represented in the form 

P{x,t) = Jl^Sit-s) 1^,-1 (\/(E)P),| ds 

= Pi{x,t) + P2{x,t) 

Recalling (p2D, we see that any bound for || A ||oo would yield at once a bound for 
II Pi{x,t) lloo- On the other hand, if || P is bounded, then {VP)x G W'^'""^ for ai 3> 1, 
whence | /q S{t — s){VP)x ds\ < Ct^ || VP ||oo for some positive C and u. A fixed point The- 
orem as that sketched to obtain local existence will then give ([l9| ) for < t < 5i <^ 1. On 
iterating in time the previous argument, (|TP|) now follows for the general case thus proving 
Lemma 2. 

We next observe that, by Lemma 1 and ( [T^ ) 

1^^^^^ l^sup |E(x,t)|"^j d?j < C {1 + T). (23) 

where C = C(||A|U,PIU). 

Taking into account the growth condition (0) it then turns out that \f{A,E)\ < 
C{\E\^ + 1) for some C = C{\\ A ||oo) > 0. By (^, we now see that 

rT+l 



1/13- 



i^J^ \f{A,E)f^dxdtj <C{1 + T), (24) 

where /?i > 1 and C = C(/3i) > 0. 

Putting together (|2^) and (|19D , we see that a = cr(ai,/3i) > 1 exists such that 
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{It^^ Io l/(^'^)^r^^^^) < C (1 + T)^ (25) 
where 6 = 6{a) > 1. In view of (|T]), At also satisfies (p3D, and this in turn implies that 

\P{x,t)\ < C{1 + tf' for some Oi > 0, (26) 

whereas 

\A{x,t)\ < C (1 + tf' for some 02 > 0, (27) 
both bounds being valid as long as P and A exist. This gives at once global existence. 

III. ASYMPTOTIC ANALYSIS 

It is convenient to redefine the time and space scales in such a way that the semiconductor 
length becomes 1: 

I X t 



Then Eqs. (i)-(|10D become 

J-J{E) 



1 



K{E) + R{E) 
din J 



-ViE) 



ds 



mi ^ (29) 

J [K{E) + R{E)] dyds ' ^ ' 

'Eiy,s)dy = ^, (30) 

E{0,s) = pJ{s), (31) 

in terms of the 'slow' variables s = et and y = ex. We shall describe the time-periodic 
solution of these equations in the limit e — (L ^ oo) by leading order matched asymptotic 
expansions. 

A. Outer solutions 

Clearly the leading order of the outer expansion of the solutions to Eqs. (^) - (^) yields 

J-J{E) = 0. (32) 

Then the outer electric field is a piecewise constant function whose profile is a succession of 
zeros of J — j{E), Ek{J) {k = 1,2,3) separated by discontinuities. These discontinuities are 
traveling wavefronts whose velocities can be found by studying the inner solutions. Typically 
pairs of these wavefronts bounding a region where the electric field is E^lJ) constitute 
solitary waves, which are the basis of our study. See below. 
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B. Inner solutions 



Near y = or near the moving discontinuities in the outer field profile, there are regions 
of fast variations of the electric field. In them J = J{s), E ~ F{x,t] s), where x = y/e and 
t = s/e according to Eq. (^81). 

1. Injecting layer 

The field in the injecting layer (defined as the boundary layer near y = 0) obeys the 
equations: 

V{F)^ / d^F JV'{F) dF\ 

J{s) [K{F) + R{F)] \dxdt ^ V{Ff 'dt j 

+ViF)^ = J{s)-jiF) 

I , ^(^) ^ (33) 
^'k{F) + R{F) ds ' ^^^^ 

F{0,t;s)=pJ{s). (34) 

Except in very short time intervals where a new wave is being formed, F = F{x; s) is 
a quasistationary monotonically decreasing profile joining F = pJ{s) and F = Ei[J{s)], 
which solves the equation: 

V{F)^ = J{s)-j{F). (35) 



2. Traveling wavefronts 

Near a discontinuity of the outer electric field profile, located at e.g. y = Y{s), the 
electric field is a traveling wavefront solution of (p3D, F = F(x;s), with x = [v ~ ^("S)]/^ 
and dY/ds = c: 



j..KiFl+RiF)_ 



cV'{F) 



K{F) + R{F) 



-ViF)]^ + Jis)-jiF)^=0. (36) 

These wavefronts should join the two different stable zeros of j{E) — J, Ei{J) and E^lJ). 
To find them we have to (i) fix J{s) = J, (ii) vary c( J) until we find a solution of (^6|) which 
is equal to Ei{J) [resp. E2,{J)] when x = "Oo and to E^{J) [resp. Ei{J)] when x = +oo. 
Thus there are two types of wavefronts, each with a given wavespeed c , which is a unique 
function of J{s): 
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A heteroclinic solution of the phase plane corresponding to (|36D joining the saddles 
(i?i(J),0) and {E3{J),0) with dF/dx > 0. For each J G {jm,jM) there is one such 
solution with wavespeed c+(J). 



• A heteroclinic solution of the phase plane corresponding to (pq) joining the saddles 
(£^3(7), 0) and (Si(J),0) with dF/dx < 0. For each J G {inn 3m) there is one such 
solution with wavespeed c_(J). 

The functions c+(J) and c_(J) are depicted in Fig. 0. Notice that they cross at a value 



J = J* which was already defined in pO] as the value of the current for which a homoclinic 



orbit joining (-Ei,0) to itself 'collides' with the other saddle (£^3,0). 



C. Putting the pieces together 

We shall start at a time when there is only one solitary wave in the sample, far from its 
ends. We shall assume that Em < (p < in Fig. |I| and that Em/Jm < P < Em/jm- Then 
a Gunn effect mediated by solitary waves occurs PT|. The initial outer field profile is 



E{y, 0) = E,{J) + [Es{J) - E,{J)] [e{Y^ - y) 

-e{Y,-y)]+0{e), (37) 

E,{J)-E,[J)- 

Here Yi and 1^2, ^1 < ^2 are the positions of the fronts connecting Ei{J) and E^{J). Fixing 
the values of (j) and J(0), (|38D fixes 12(0) — ^1(0). 9{x) is the unit step function: 9{x) = 1 
if X > 0, 9{x) = if X < 0. The outer field profile for s > is given by (|37| ) with J = J{s), 
and Yi = Yi{s) {i = 1,2). The location of the wavefronts is given by the equations: 

whereas their separation is related to the bias through (BSf). We can find an equation for 



J by differentiating ( p8|) and then inserting (pOf ) in the result. We obtain: 

dJ (E3 - E,)' 



ds + 



(c+-c_) (40) 



where j- = j'[Ei{J)]. This is a simple closed equation for J that demonstrates that J tends 
to J* (for which c+ = c_) exponentially fast. After a certain time, the wavefront Y2 reaches 
1 and we have a new stage governed by ( |37D with Y2 = I and Yi given by (j3^). The equation 
for J becomes 

dJ (Es-Ei)^ , , 

and J increases until it surpasses the value such that E^i^J) = pJ at some time s = Si. 
After this time the boundary layer near ?/ = is no longer quasistationary. In fact, about this 
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time a new stage begins during which the boundary layer becomes unstable (see Section |VD 
and it starts shedding a new solitary wave on the fast time scale t = {s — Si)/e. Meanwhile 
the current density reaches its maximum and it decreases until J = Jc again later on. During 
this stage the new solitary wave is generated and it starts departing from y = 0, a process 
which will be described in the next Section. Let us call T2 > the fast time at which the 
solitary wave is fully formed, so that the field inside it is E^{J), and outside it is Ei{J). The 
wavefronts enclosing the high field region are centered ai y = and y = Y4, Y3 < Y^. We 
then have another slow stage, during which the outer electric field profile is: 

= E,{J) + [Es{J) - E,{J)] {1-Y, + Y,- Ys) 

+0(e), (42) 



(assuming that (f) is so large that the old wave Yi{s) < 1). Differentiating (|42|) with respect 
to s and using that Yi and Ys move with speed c+(J) while I2 moves with speed c_(J), we 
obtain: 

dJ (E^-E^ 

^ = ^.E, I Es-^ (2c+ - c_) . (43) 

We have now two possibilities: either (i) (2c+ — c_) < and then J further decreases until 
the zero of [2c+(J) — c_(J)] is reached (again for biases so large that the old wave at 
Yi{s) < 1); see Fig. H, or (ii) (2c+ — c_) > and new wave(s) are shed; see Fig. |[ We 
consider only the stable case (i) which occurs for Jc > large enough (or equivalently, for 
small enough contact resistivity p). In this case J{s) decreases exponentially fast towards 
J^, thereby forming a plateau when the bias (j) is large enough. After the old wave reaches 
y = 1, we get again Equations (|37| ) - (|iO| ) with Y3 and Y4 playing the roles of the previous 
wavefronts Yi and Y2, respectively. We thus recover the initial situation and a full period of 
the Gunn oscillation is completed. 



IV. SHEDDING NEW WAVES FROM THE INJECTING LAYER 

After the current density surpasses the value Jc at s = Si, a new stage begins during 
which a new wave is shed from the injecting layer. To understand how the quasistationary 
injecting layer solution becomes unstable and a new wave is shed, it is convenient to keep 
two terms in our asymptotic outer expansion of the solution: the second order term in the 
outer expansion enters at the same order as the leading order term of the injecting layer 
in the bias condition. It can be seen that the stability of the injecting layer is governed 
by equations similar to those of the stationary solutions |T^J2B| , except that the old wave 



located ai y = Y{si) needs to be taken into account in the outer quasistationary profile. 
Then a similar instability criterion holds, and the injecting layer becomes unstable when 
F ~ E2{Jc) on an interval of width O(elne) next to y = 0. 

The appropriate time scale for the shedding stage is r = (s — Si)/e, so that position of 
the old wavefront is Yi ~ Yi(si) + ec+r during this stage. Inserting the Ansatze 

E = E('\y,T) + eE('\y,T) + 0{e'), 
A = A^'\y,T) + eA^'\y,r) + Oie'), 

j = j^ + eJ^'\r) + 0{e^), (44) 
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in Eqs. (|^) and (H), we find that the leading order term E^^\y,T) is given by (|37|) with 
Yi ~ ^i(si) + ec+r and F2 = 1, and J = J^. After some algebra, the equations for E^^^ and 
Ad) turn out to be 

i?« = ;^(J«-V^.^«), (45) 

= ^i^a'V^^-^cV/)^^^^- (46) 

In these equations the functions of the outer electric field profile are piecewise constant, for 
example, V = V{Ei) = Vi if < y < Yi{si) and V = 1/(^3) = I/3 if Yi{si) < y < I. The 
subscripts i are either 1 or 3 according to the value of y. The solutions that do not increase 
exponentially as r ^ —00 (and thus they may match the previous stage) are 



XI e ^'J^^\t)dt. (47) 

-00 



Here 

Iii = , « = 1,3, (48) 

i = liiO<y< Yi(si) and i = 3 if Yi(si) < y < 1). The expression for A^) immediately 
follows from (^) and (|7). 

The equations at the injecting layer are (0) and (|) with J = + e p J^^\t) + O(e^), 
while the boundary condition is F(0, r) = E2{Jc) + ep J^^K^) + 0{e'^). To determine J^^\t) 
we shall substitute the outer field profile ( ^41) in the bias condition (H). The result is 

= Ei(Je) + [E3(^c) - ^l(Jc)] [1 - Y,{S,)] 

+6 {4'^ + [E« - Ei^)] [1 - Y^{S,)] - [E,{J,) 

-E,{J,)]c+T + / [F{x,r) - E,{J,)] dx 
Jo 

+ f [F{x)-Ei{J,)]dx 
J —00 

+ flFix) - EM] dx] + 0{e'). (49) 



Let us define now 



PCX) 

h{T) = [EsiJc) - E,{J,)] c+T - / [F(x, r) 

Jo 

-Er{J,)]dx- f [F{x)-E,{J,)]dx 

J —00 



I —00 

coo 



POO 

/ [F{x)-E,{j;)]dx, (50) 
Jo 
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which is essentially the area lost by the motion of the old wavefront Yi during the time r 
minus the excess area under the injecting layer. If we include in the leading order of the 
outer electric field the two terms of (|i9| ) that represent the excess area in the transition layer 
at Yi, and insert (^Tj) in the result, we obtain 



^1^1(51) ^J/3[l-n(Sl)] 



JcVi 

-Jcv; 



JcVi 



^3 [l-Yi{si)] 



J —00 



(51) 



After some elementary manipulations, we can transform this integral equation for J^^^ 
in a linear second order ordinary differential equation which can be readily solved under the 
condition that the solution should not increase as r ^ — 00. We obtain: 



T 



ai 



G{a) ^^"-"^ da. 
ai 



(52) 



^ (Ai + /ii)(Ai + /i3)e^^" - (A2 + /ii)(A2 + /i3)e^^" 

""^"^ = aT^a^ • 

Here Ai and A2 are the (negative) roots of the quadratic polynomial: 

OiA^ + a2A + 03, 



whose coefficients are 



_ ViY.jsi) Vs[l-Y,{s,)] 

JcV{ + J.Vi ' 

(Y^{s,) V,[l-Y,{s,)]\ 
— f'l I r 



J'l 



+/W3 



JcVi 

yiri(si) , i-Fi(si) 



J3 



JcVi 

(Y^{s{) , i-ri(si; 



3[ 



(53) 



(54) 



(55) 



(56) 
(57) 



We can write (|52D in another form that suggests a more transparent interpretation: 

(r) = J'(si) (r - ro) - /(r) - / G(a)/(r - a)da 

Jo 

= /(si)(r-ro)-/(r), 
where J'{si) is given by (|l|) at J = Jc, G{a) is the kernel (|53|), and 



(58) 
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/(r) 



-1 



1 



To 



(E3 - Ei)c+ 



+ 



[F(x)-i?3(^c)]rfx 
il 



[F(x,r)-Ei(J,)]rfx, 

r [F(x)-i5;i(Jc)]rfx 
J —00 



(59) 



X 



(E3 - Ei)c+ 
Yr{s^){j[Vi-J,Viy 



JcfJ'lVlj'i 



(60) 



The terms on the right side of (|58D clearly display the balance between the area lost by the 
motion of the old wavefront at ^1(5) and the excess area under the injecting layer. 

Let us recapitulate. We need to solve the following equation for the field in the injecting 
layer: 

d'F ^ [J, + eJ ^'^T)][K{F) + R{F)] 



dxdr 



V{Ef 

K{F) + R{F) dt ^ ^ dx 



or equivalently, the system 

dA 
a7 



V{F) dr 



[Je + ejW(r)] [KiE) + RiE)] 
V{Ey 

x[j{E)-ViE)A], 
_ J^ + eJ^^\T 

dx 



(61) 



V{E) 



A, 



for A and E = F{x, r), with the following boundary condition 

F{0,T) = p[J, + eJ^'\T)]. 



(62) 
(63) 

(64) 



Here J^^\t) is given by (|52[) as a functional of the function /i(r) or by (|58D as a functional of 
the excess area under the injecting layer, /(r). Besides solving this problem, F{x,t) needs 
to satisfy the matching condition: 



F(a;,r)-Fo(a;;J(s))<e, 



(65) 



as s — i> si — , r — > —00 on an appropriate overlap domain. Here Fo{x; J{s)) is the quasis- 
tationary boundary layer solution of (|35D and (|3^) for s < si, J{si) = Jc, with J(s) given 



by (prp. We have kept the 0(e) term in ( p^ because it becomes of order 1 after the new 
wavefront (s) are formed. Eqs. (|6TD - (^) constitute a meaningful problem which can be 
solved numerically or further studied by analytic means. Notice that the forcing term at the 



13 



boundary, proportional to J^^^(r) is a functional of h^r): the area lost as the old wavefront 
moves toward y = 1 minus the area gained because of the field growth in the injecting 
layer. The present discussion could be used to make precise the observations of Higuera and 
Bonilla for the creation of a new wave in the usual model of the Gunn effect (see Section 5 



of [24|). Anticipating the creation of new wavefronts near y = 0, we decide to define their 



position on the slow space scale by Yi = eXi {i = 3, 4), X3 < X4, such that 

f(.Y3(r),r)=e„:=:?ji:^^l±^, (66) 
F(Jf4(r),r) = £„:=M^=l±^. (67) 

The results of our simulation are shown in Fig. ^ When -F(x, r) is sufficiently close to 
E^[J{t)] on some space interval, we consider that the fast shedding stage has ended. Then 
the positions of the new wavefronts ^3 and 14 and the value of the current density may be 
used as initial conditions for the following slow stage with three wavefronts as described in 
the previous Section. 



A. A local limiting case when In L » 1 



There is an asymptotic limit in which Eqs. (|6l|) - (|65|) can be recast as a simple, local, 
universal problem. On the r time scale and whenever J^^)(r) < 0, the quasistationary 
boundary layer solution of (^) and 
satisfies (see Appendix A of [pH]): 



[as (x — X) 



and 



[as (x — X) 



Fq{x]J{s)) E2- cl exp 

-00], with 

cl = {E2 - Eo) exp 

J2 V{E) 



for s < si, J{si) = Jc, with J(s) given by (^if) 

J2 



Vo 



E2 
Eo 



V2[jiE)-J] 



E — E2 
dE\, 



Fo{x; J{s)) Ei + cr exp 
-ool, with 



J'l jx-X) 



Eo 



Cr = {Eq - El) exp I - ^ 
fi V{E) 



Vi m) - J] 



E-Ei 
dE\ 



(68) 



(69) 



(70) 



(71) 
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Here we have defined X4 = X{s) such that Fq{X{s); J{s)) = Eq given by (|67|). By using the 
boundary condition at x = 0, (^) and the matching condition (^), together with (|6^) , we 
obtain 



•'2 



X 



J'{s,){t-to)-I{t) ^-CLe—, (72) 



which yields the instantaneous "position", X{s), of the transition layer linking E2 and Ei 
as s — > — Si. Of course during the shedding stage the position of the transition layer is no 



longer given by (|7^). To find it, we need a better approximation to the field at the injecting 



contact than the quasistationary one. We shall proceed as follows: 

• We define (as before) X{t) as the position for which the field in the injecting layer is 
(El + E2)/2. 



We find an approximate solution to Eqs. (|6T|) - (^) in the limit Ax ^ — 00 where 



X = X(r) + Ax and 1 <^ |Ax| <^ X(t). This approximation breaks down if dX/dr = 
0(1), which occurs for r = O(lne). This growth of dX/dr indicates that the end of 
the injecting layer starts advancing and shedding of a new pulse starts. 

For those times, we can find a reduced shedding problem with a local boundary con- 
dition. Its solution describes the formation of a wavefront joining E^{J^ and -Ei(Jc) 
(which will eventually become the foremost wavefront of the newly shed pulse, ^4(5), 
described in Section |T|). After a long time, the back of the pulse is a stationary layer 
joining E2{Jc) and E^{Jc). 

The formation of the back of the new pulse will be described by an equation similar 
to the previous reduced shedding problem with different local boundary and matching 
conditions. 



1. Approximate form of the electric field near the injecting contact: forefront shedding 



Far from x = X(r) and close to x = 0, we can linearize Eq. (|6lD or equivalently, Eqs. (^ 
and (|5|) about A = JjV{E2{Jc)) and E = EiiJc). The result is 

d'^E [K2 + R2] dE J^dE_ 
dxdr V2 dx V2 dr 



Jcf2 [K2 + R2 



1 djW 



'' = '\V2 dr 
E(0,r) = epjW(r). 



(73) 
(74) 



We have ignored terms of the second order in e and the linearized field. The solution of 
Eqs. (|7|)-(|7|) which satisfies 



Eix, t) e 



' \\ _£k 1' 
P--7 e +- 

\ J2 I J2 



J«( 



T 



(75) 
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(as r — oo) is derived in Appendix 0. Its asymptotic behavior as x — > +00 corresponds 
to the transition region where x ~ X{t) ^ 1: 

E^e(^p-j^ e"^ {J'(si) [r-To-ex] 

-i{r-ex)}, (76) 

B - ^'^^^-^'^^^ > n (77) 

To discuss the hmit of vahdity of (|76D , we consider the region where E = 0(1). For instance 
we may set 

e(p-\] e-'-W [J\s,) (r - 9X) - /(r - 9X)] 



J2, 



(78) 



according to the definition of X{t) as the position for which the electric field takes on the 
value {El + E2)/2. When r < and |r| ^ x and using the fact that | lne| 3> 1, we just 
neglect the second term on the left side in ([78|) to obtain 

X(r) ~ ii m (^El^3}2L\ . (79) 

This approximation (the adiabatic approximation) holds if dX/dr <^ 1. Our approximation 
continues to hold even when r increases to zero in such a way that |r| ^ X. Actually the 
adiabatic approximation breaks down when the left side of (|7^) changes sign. Once this 
happens, we can no longer discard the term proportional to /(r — 9x) in ([78|) . We have 

i{T -6x) ^ K J\si)X{t ~eX{T)), (80) 

K is approximately constant. Then the breakdown condition is 

r - eX{r) ~ r^KXir- 9X{t)), (82) 



whence 



.2(l-pj0W(si)e 
We now set new length and time scales as follows: 
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(85) 



IJ2I V2(l - PJ2)V2eJ'{sr)e^ 

While this change of coordinates leaves invariant (|^) and (|^), inserting condition (^) in 
(^), we find 



x{a — 9r])e^, as 77 — — 00, (86) 



We are thus led to the following reduced shedding problem. Solve for — 00 < < cxd, 
and —00 < (T < 00: 

(87) 

subject to the boundary condition ( pB]) and to the matching condition ( |65D (once the changes 
of coordinates (|8^ -(p5D are inserted in it). 

The boundary condition (|86|) is local and it represents a small wave whose amplitude 
increases exponentially as it moves to the right with velocity 1/9. Notice that the same 
interpretation can be made of the solution of the amplitude equation describing the onset 
of the Gunn-type oscillations in the supercritical case [|1^. Unlike the bifurcation case, the 
solution of the problem (|86| ) - (|88|) describes how the amplitude of the electric field wavefront 
grows until it reaches E = E^lJc), which happens at a time we will denote by r*(e). Then the 
injecting layer presents the following structure: a forefront which is an advancing transition 
layer joining E^[Jc) and £^i(Jc) and located at X(r), a region of width (Ine) on which 
E = Esi^Jc), and a quasistationary back front which is a transition layer joining E2{Jc) and 
E-si-Jc). 



2. Approximate form of the electric field near the injecting contact: backfront shedding 

Clearly what happens next, for r > r*, is described by other equations. Here we shall 
only consider the stable case, Jc > J\ in which only one wave is shed. In particular the 
expression for /(r) is different and J^^^t) starts decreasing so as to keep the bias condition. 
There is another local problem which governs the transformation of the quasistationary 
backfront joining E2{Jc) and E^lJc) into a moving wavefront joining -E'i(Jc) and E^{Jc). 
Since the necessary analysis is similar to that for the forefront, we will describe only its 
more salient features, omitting many details. 

Let us denote by X^{t) the instantaneous position of the forefront from now on, and 
let X{e) = Xi{T*). Let us make a new change of coordinates which is simply a shift of 
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X = X{e) + r/, (89) 
r = r*(e)+(T. (90) 

The excess area in the injecting boundary layer is no longer proportional to {E2 — E1) X^. 
We have instead that for r > r*, this area is equal to the sum of: (i) the excess area at the 
time r*, (ii) the excess area due to the motion of the forefront between the times r* and r, 
minus the loss of area due to the motion of the backfront located at X^la): 

a,Iia)= / [F{x,cr)-E,{J,)]dx = iE2~E,)X 
Jo 

+ {E, - E,) c^a- {E, - E2) [X - X,{a)] (91) 

X3(cr) is defined as the position for which the electric field takes on the value {E^ + E2)/2. 
Then 

/ 1 \ ^2-^3 

e e~—[J'{s,){T-eX,) 



J2, 



E^ — E2 



-/(r* + a - 9X,)] ~ (92) 
and therefore we can make an approximation similar to ([791), 

A-3M~||ln(^)~A--i|l„., (93) 

for times a < 9X. F > is a constant which we shall leave unspecified. After times a = 6X, 
the approximation given by (0) breaks down and / in (^) should be replaced by (0). The 
adiabatic approximation breaks down when the left side of (^) changes sign, which yields 
new times a* and positions X3, after which the backfront is shed. To analyze the shedding 
problem we define new variables a and rj by using ( |5UD and (|^) with a* and X3 instead of 
T* and X. The result is that the shedding of the backfront is governed by Eqs. {^7\j - (|88D 
with the following boundary condition instead of (^): 

E{r], a) ~ E2{J,) - _ 0^) (94) 

Zt 

(as T] —00; r' is a positive constant), and the new matching condition with the core- 
sponding stationary state. 



V. CONCLUSIONS 

We have shown that a unipolar drift- diffusion model of charge transport in ultrapure 
p-Ge has solutions existing globally in time. We have performed a new asymptotic analysis 
of a reduced version of this model which describes the trap-dominated slow Gunn effect 
on a long sample. The building blocks of this analysis are the heteroclinic orbits used to 
construct the typical solitary waves mediating Gunn-like oscillations. During most of the 
oscillation, the motion of the heteroclinic orbits and the change of the electric field inside and 
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outside the solitary waves (enclosed by heteroclinic orbits) follow adiabatically the evolution 
of the total current density. When a solitary wave reaches the receiving contact, the current 
increases until the injecting contact layer becomes unstable and sheds new wave(s). The 
wave shedding stage of the oscillation is described by a semi-infinite problem which needs 
to be solved and matched to the rest of the oscillation. As an outcome, we have found a 
criterion that shows that single or multiple wave shedding is possible during each oscillation, 
depending on the resistivity of the injecting contact. While single shedding is the usual 
(stable) Gunn effect, multiple wave shedding breaks the spatial coherence of the sample. 
This new instability mechanism may eventually explain the complicated behavior observed 
in experiments performed in long semiconductor samples We have confirmed these 

results by direct numerical simulation of the reduced model, in particular the new predictions 
of multiple shedding of solitary waves in the unstable case. A simpler asymptotic description 
is possible in the limit InL — oo: we need to solve two spatially infinite problems, one for 
the shedding of the forefront of the new solitary wave, the other for its backfront. Although 
this limit is unrealistic and extremely hard to compare with numerical simulations {L should 
be so large that the needed computing time is excessive), we have analyzed it in some detail 
in view of the insight offered in the shedding problem. 

We expect that a similar analysis of the Gunn effect can be performed in other models, 
although the details of the shedding problem will in general be different . 
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APPENDIX A: SOLUTION OF THE PROBLEM (73)-(75) 

In the sequel we shall study the form and asymptotic properties of the solutions of 
i|)- ([75|) . Eq. ( [73D may be written as 

OXOT OX OT 



V2 [ dr 



(Al) 



J^V^ J^f, [K2 + R2] 



^ ^ MK2 + R2) ^ Bf,2V2 _ ^^2) 

^2 j2 



We now write E{x, r) as a sum of the solutions of simpler problems: 
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i=l 



CSi = 0, £i(0,r) =pJ'(si)(r-ro), 

J'{s,)[l + Dir-ro)] 

, ^2(0,r) = 0, 



Vo 



1. 



C£3 = 0, ^3(0,r) = -p/(r), 
+ Di] -gir), A(0,r) = 0. 



(A3) 
(A4) 
(A5) 
(A6) 
(A7) 

(A8) 
(A9) 



where 9 is given by ([77|). To calculate £3, we observe that, on taking Fourier transform in 
f), we obtain 



The solutions of ( [A^ ) and ( |A5|) are immediate: 



£1 = p J'(si) (r - To - e x) e , 

DJ'{si) ^ . Qx 

^2 = [(r - To - ^x) e « - (r - Tq)], 



£^3(x, r) = — p / K{x, T — s) I{s) ds, 



where 



K{x, r) = — K{x, Lu) e"^^ du, 

ZTT 



exp 



{DB + C)x' 
D + iuj 



- 1 



(AlO) 

(All) 
(A12) 



Hence 



r 1 

K{x,r) = e-^- [5{r) + - j 

-1] e'^^ duj] = h{x,T) + h{x,T) 



{DB+C)x 



(A13) 



The term /2(x,r) in ( |A13|) can be evaluated by standard complex variable methods. To 
this end, we notice that the residue of the function expi(Cr — M/C) at C = is 

1 g M'+V' 



Then Cauchy's theorem yields 

J2(x,r) = e"^"-^" (D5 + C)a; 

x^[(D5 + C)xr]77(r), (A14) 

where //(r) = 1 for r > and //(r) = otherwise. The function ^ is defined as follows 

_ /i(2Vi) 



;e^/!(/ + l)! 



(A15) 
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where Ii{z) is the standard modified Bessel function of the first kind with index 1. From 
(|AT0|)-(|M5D, we obtain 

£s{x,t) = -pe-^" {iir) + {DB + C)x 
X r e-^'^^-'^^[iDB + C)x{T~s)]i{s)ds]. (A16) 

We now turn our attention to S4. We claim that 

S4{x, t) = F{t) - r K{x, T-s) F{s) ds, (A17) 

J —00 

where F{t) is any arbitrary solution of the equation 

dF 

B — -CF = g{T), (A18) 

with g{T) as in ( |A7|) . At this juncture the reader may wonder whether (i) S4 thus defined 
is actually independent of the choice of F, and (ii) £4 is such that £4(0:, r) depends only 
on those values of g{s) for which s < r. It turns out that the answer to these questions is 
yes. However, not to interrupt the main fiow of the arguments here, we shall postpone the 
discussion of this point till the last part of this Appendix. We continue by describing in 
some detail the asymptotics of £^3(0:, r) and £^4(x, r) as x ^ 00. This is all we need to get 
the asymptotic form of E in view of the explicit form of Si and S2 given by (^5) and ( [A 51 ) , 
respectively. 

Let us begin by S3. Since by ( [77D and ( |A^ ) 9 = {B + C/D)/D, we need to study the 
behavior of the second term in the right side of ( |A16|) : 

R{x, t) = -pD^Oxe-^'' r e~^("-^) 

xm[D'^ex{T - s)\i{s)ds, (A19) 

as X — > 00. Without loss of generality, we may assume that x{t — s) ^ 1 and make use of 
the approximation Ii{z) ~ (27r2:)~^/^e^, as 2; — 00. Then 



pOxe "-^ f 



2v^ ^-oo [D^dx{T-s)]i 



™-(2v/i-0rf^. (A20) 



xe 

In view of its definition (|^) and (0), it is natural to assume that /(r) has algebraic 
behavior as r — > 00. Then the main contribution to R{x,t) will correspond to ^ ~ 1, and 
we may approximate the integral by Laplace's method, thereby obtaining 

S-i{x,T) R{x,r) pe^i{T-0x) (A21) 

(as X ^ 1), after substituting the result in ([A16| ). Let us consider now £^4(0;, r) given by 
( |A17|) . By selecting the function F having a zero at r, F{t) = 0, we observe that 
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D 



CVo 



K{x, T — s) F{s) ds 



e o I{t — Ox), as x ^ 1. 



(A22) 



Putting together (^), (|A|), ([MID and (|A22D , we finally get 



-i{T - Ox) 

as X ^ 1. It is easy to check that this equation is the same as 



(A23) 

by using the definitions 



We conclude this Appendix by stressing that S^i^x, r) given by (|A17|) is indeed well 
defined. This will be consequence of the following Liouville-type result, which is of some 
independent interest: Let u{x, r) be a solution of 



— — + D — + B—-Cu = 0, 

UXOT OX OT 

where D, B and C are positive constants, and assume that u is such that 

n(0,r) = 0, 

\u{-,t)\ < u[x) e'^^'^\ as r ^ oo with a < D. 



(A24) 



(A25) 



Then u{x,t) = 0. To check this, we set r = — s in ( |A24|) and take the Laplace transform of 
this equation, thereby obtaining 



iD-z)^-{Bz + C)U 



du{x, 0) 

dx 
-Bu{x,0), 



(A26) 



for U{x, z), the Laplace transform of u{x, r). Assuming that (|A25| ) holds, U will be analytic 
in Im^; > 0, and integration of (|A26|) yields 



U{x,z) 



1 



D — z Jo 



Bz + C 



X 



dx 



D-z 
-Bu{^,0)] 



(A27) 



Then U{x,z) as \z\ oo, and therefore U{x,z) = 0, whence the result. To derive the 
desired result for S^^x, r), we merely observe that 



f K{x,T - s) es" ds = eB'^. 

J —oo 



(A28) 



Indeed the function on the left-hand side above is a solution of ( [A24| ) , and so is es'^. Since 
the difference between these two solutions is at most algebraic as r ^ — oo, ( [A25| ) holds and 
the result follows. 
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Using the fact that two different solutions of (|A18|) differ by a term ki e^'^^^ (for some 



real constant /ci), it follows from ( |A17| ) and ( |A28| ) that £4 is independent of the particular 
choice of F. On the other hand, if we select F (for a given r = f) as follows 



F{s;f)= I ei(^-«(7(0cie, (A29) 



then it turns out that for s < t, F only depends on the values of /(^) with ^ < r. By (|A17|) 
the second claim on S^i^x, r) is thus proved. 
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FIGURES 



FIG. 1. (a) Drift velocity, (b) impact ionization coefficient, (c) recombination coefficient, and 
(d) homogeneous stationary current density j{E) for a = 1.21, illustrating the NDR in our model. 

FIG. 2. Velocities c+ and c_ of the heteroclinic wavefronts as functions of the current density. 
Notice that the lines intersect at J = J* ~ 0.076515. We have also marked the current ~ 0.07840 
at which 2c+ = c_. ji^^ and ji^^ are the critical currents Jc for the resistivities of Figures 3 and 
4, respectively. The inset shows the location of J^, ji^^ and ji^^ on the bulk current density curve 

m- 

FIG. 3. Numerical solution of the reduced model during one period of the oscillation in the 
stable case 2c+(Jc) < c_(Jc) for which only one new wave per period is shed. Parameter values 
are L = 3000, p = 7.49736, = 0.749565. 

FIG. 4. Numerical solution of the reduced model during one period of the oscillation in the 
unstable case 2c_)_(Jc) > c_(Jc), corresponding to p = 10 and the rest of parameters as in Figure 
3. With these values of the parameters, two new waves are shed during each period. 

FIG. 5. Comparison between the results of direct numerical simulation and the leading-order 
asymptotic approximation of the solution for the same values of the parameters as in Figure 3. 
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FIG. 2. 
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